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ABSTRACT 



Heat transfer rates for laminar, convective heat transfer in the 
entrance region between parallel plates are investigated. The hydro- 
dynamic solution due to Bodia /2/ was used in the solution of the energy 
equation in finite difference form on a digital computer. The thermal 
boundary conditions include: constant heat input, constant wall tempera- 
ture, one wall constant temperature and one wall insulated, and constant 
but different wall temperatures on the upper and lower walls. 

The approximate, integral methods of Siegel and Sparrow/^, /?_/ 
produce results that are in close agreement with the solutions in this 
analysis for the constant heat input and constant wall temperature cases. 

The scope of the finite difference solution is limited to a narrow 
range of Prandtl numbers near unity, due to the small grid sizes required 
for convergence at small Prandtl numbers and to the overly low transfer 
rates indicated near the entrance for high Prandtl numbers, which is a 



result of the "finite starting length," 
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TABLE OF SYMBOLS 



English Letter Symbols 

Specific heat 



c 

d 

D 

h 

k 

L 



■•m 



w 



W]^ 



W2 



T 

T 

I 

u 



w 



m 



Half the distance of separation of the parallel flat 
plates 

Hydraulic diameter D = 4d 

Unit conductance for convection heat transfer 



Unit thermal conductivity 
Dimensionless length parameter, 
Jleatrflux'.perl.unit area.; qi/A 

qD 



7d = 



Re, 



Akt, 



X 



16 



Dimensionless heat flux, 

Temperature 
Initial fluid temperature 
Mixed mean temperature 
Wall temperature 

Lower wall temperature (asymmetric case) 

Upper wall temperature (asymmetric case) 

Dimensionless temperature, t/t^ 

Dimensionless wall temperature, t^/tj 

Dimensionless mixed mean temperature, t /t^ 

m o 

Axial velocity component 



X 





Initial axial velocity component 




U 


Dimensionless axial velocity, u/Uq 




V 


Crosswise velocity component 




V 

c 


Crosswise velocity component at centerline 






Initial crosswise velocity component 




V 


Dimensionless crosswise velocity component, v/uq 


X 


Dimensionless length parameter, 

d^Uo 




y 


Normal coordinate 




Y 


Dimensionless normal coordinate, y/d 




Greek Letter 


Symbols 






Fluid density 






Fluid viscosity 




n) 


Kinematic viscosity, 




Non-dimensional Groups 






Local Nusselt number, 




PR 


Prandtl number, 

k 




®d 


Reynolds number, based on half the gap width. 


dUo 


V 


®D 


Reynolds number, based on hydraulic diamjeter, ■ 


Duo 


V 



xi 



i 






f 

-If 






I. INTRODUCTION 



Efficient design of compact heat exchangers such as are employed 
for gas turbine regenerators, among other applications, requires maximum 
utilization of the high transfer rates available in the entrance region. In 
this region the velocity and temperature profiles of the entering fluid are 
undergoing rapid transition from their uniform distribution at entrance to the 
fully esta()lished distributions encountered farther downstream. The well 
known Graetz solution, in which only the temperature development is con- 
sidered, does not yield adequate results when the velocity and temperature 
distributions are simultaneously developing. This is particularly true if the 
entrance region represents a considerable percentage of overall length and 
the two profiles develop at approximately the same rates „ This situation is 
frequently encountered with gases where the Prandtl number is approximately 
one. 

The purpose of the present analysis is to present a finite difference 
solution for the heat transfer rates between parallel flat plates for several 
different boundary conditions corresponding to various constant wall 
temperatures and constant heat input configurations. The approach to the 
problem was suggested by Miller l_\J , which is adapted from the numerical 
procedures applied to the entrance region of the circular tubes by Kays 1_%J . 
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While an analysis such as the present one can hardly compete for 
simplicity with approximate treatments such as the Graetz analysis, or 
the several energy integral analyses commonly employed in engineering 
design, it does provide an exact solution for a limited number of cases 
which serves as a standard with which the accuracy of the various approxi- 
mate methods may be assessed. 
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II. ANALYSIS 



A. Governing Equations 

The governing equations are the energy, momentum (Navier-Stokes) 
and continuity equations . In order to reduce the complexity and coupling 
of the equations, the following assumptions are made. 

The flow is assumed to be: 

1 o steady, two-dimensional 

2 . laminar 

3. incompressible 

In addition it is assumed that: 

1. thermal diffusivity ) is constant 

2. convective heat transfer is large compared to radiation, axial 



conduction, and viscous dissipation. 

The governing equations in reduced form may then be expressed as: 
Energy: 








( 1 ) 



Momentum: 




( 2 ) 



Continuity: 






(3) 
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Due to the assumption of incompressibility the above equations are 
no longer coupled. The temperature development remains dependent upon 
the velocity profile , although dependence of velocity on temperature has 
been removed . 

The boundary conditions associated with equations (l) , (2), and 
(3) for flow between parallel flat plates with uniform velocity and tempera- 
ture at entrance are: 

1 . At the entrance ( x = 0) 




2 . At the walls ( y = + d) 

V = v^ = 0 
u = u^ = 0 
t = t^ 

3. On centerline 

V = Vc = 0 

4. Applied thermal boundary conditions 

(1) symmetric cases (with respect to duct centerline) 
a . constant wall temperature 



4 



b. constant heat Input 




(2) asyiBtroetric leases 

a. one wall constant temperature, one wall insulated 




b„ constant wall temperature 




where: 

Solution of the energy equation requires that a hydrodynamic 
solution be first obtained . 



B. Hydrodynamic Solution ^ 

Several hydrodynamic solutions for the entrance region are available. 
Solutions have been presented by Bodoia Schlichting /_ , Schiller /5 _/ , 

i 

and Han /6_/. The approximate solutions of Schiller and Han have the 
advantage of being analytic, and therefore, less cumbersome to use. They 
are, however, inherently less accurate than the numerical solutions presented 
by Bodoia and Schlichting. 
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n 

Bodoia's solution was obtained by evaluation of the Prandtl boundary 
layer equations in finite-difference form. Velocity profiles were obtained 
at given streamwise locations by utilizing matrix methods to simultaneously 
solve a column array of momentum expressions normal to the flow. 

A comparison of Schlichting's series solution and the Bodoia solution 
indicates that, in addition to the velocity gradient discontinuity in the 
Schlichting solution where the upstream and downstream solutions are joined, 
Schlichting's representation of the pressure gradient, based only on the 
centerline velocity, results in a core velocity gradient which is too large 
and a wall velocity gradient which is too small /3/. 

Comparisons of the dimensionless pressure decrement versus length 
due to Bodoia, Schlichting and Han are presented in Fig. 2. Velocity dis- 
tribution comparisons between Bodoia and Schlichting, and Schlichting 
and Schiller, are shown in Figs. 3 and 4. 

The Bodoia solution was adopted for the present analysis. 

C . Energy Solution 

Equations (l) and (3) may be written in dimensionless form by 
introducing the following dimensionless variables: 

T = t/t X = = x/d u = u/u_ 

o d^uo Rg ° 

Y = y/d V = v/Uq 
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In dimensionless form, the equations are then: 



Energy: 

V iT+ReV ^ = -L 

PK 






Continuity: 






+ 



Re ^ = 0 






( 4 ) 



(5) 



The following finite difference approximations are then introduced: 



- AT _ T xfijY ~ T KjV 

AX AX 



_ Tx,y»i ~ Tx.v'-i 
■ AY ■ 2 AY 



~ 2 Tx.y T Xj Y-i 

ay' 



1IL= AH = Uy.I.Y - Ux.i,Y 

AX 2 AX 



( 6 ) 



(7) 



( 8 ) 



(9) 



= AV s V )f.Y ~V x.Y-l 
AY ay 



( 10 ) 



Note that (6) is a forward difference and (10) is a backward 
difference, whereas the remaining equations are of the central difference 
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I 



I 

■ 



i 

t- 




K 

i 



type. These have been expressed in this manner for Simplicity and 
convenience. The larger truncation error normally associated with the 
forward or backward difference is not significant here since practical grid 
dimensions require that: 

AY » AX 



and, due to the nature of the flow: 

A17 \\ AV 
A X ^ Y 

Introducing (6) through (10) into (4) and (5) and solving for the 



appropriate term yields: 






AX 

trPRAY* 



+ 



AXReVl -r 
U2AYJ 







2 AX 

U PRAY^ 



] 



T 



x,Y+l 






r AX 

L UPRAY® 



AXReV 1 j 
U 2 AYJ 



«,Y+I 



(11) 



Vy,Y ■ Vx,Y-l + -R CeIx [ (12) 

which constitute the computing equations . 

Note that the Reynolds number (R^) , although appearing in the 
above equations, remains only for convenience and could be eliminated. 
It is no longer a parameter in the set of equations having been included 
in the dimensionless length variable (X) . 
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way: 



A mixed mean temperature (t^^) may be expressed in the following 



t 



m 



r cp ^ u t dA 
J" Cf, ^ u c/A 



(13) 



Removing the constants from within the integral, and simplifying, a 
dimensionless mixed mean temperature may be written in finite difference 
form: 



S Ux,yTx,y 


Z U,., 



(14) 



Due to the no-slip condition imposed upon viscous flows, heat 
transfer to or from the fluid at the wall may be expressed in terms of 
conduction through a thin film of stagnant fluid. Thus: 

f 

q" = q/A = k(||)^,„ (15) 

which must equal the heat flux due to convection, 

q" = q/A - h(t^ - t^) (16) 

Defining the local Nusselt number as; 

07 ) 

Combining (14) and (15), 
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i 

{“ 

I 

\ 

I 

I 



A Nusselt number based on the hydraulic diameter (d) instead of 



half the gap width (d) in finite difference form becomes; 



N 




(19) 



u. 



X 



T* m 



The hydraulic diameter was introduced in (19) to enable comparisons 
with previous solutions. The conversion is: 



In addition, for the asymmetric constant wall temperature case, it 
is convenient to use a dimensionless parameter other than the Nusselt 
number. A dimensionless heat flux parameter may be defined by using 
(15) above: 



D . Thermal Boundary Conditions In Finite Difference Form 

Constant Wall Temperature fe~mmetrlc Case) Since the temperature 
distribution is symmetrical with respect to the centerline, computations 
were performed in only the upper half of the duct. A coordinate system was 
adapted with the X axis on the centerline of the channel. The thermal boundary 
condition at the wall was applied as: 



D = 4d 
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1 






N 



where N is a positive integer. 



Constant Heat Input fifyintAfttijic Case) With the axis system as above , 
the slppe of the thermal profile was maintained constant at the wall by 
adding a constant to the calculated temperature one step from the wall, 
thus: 



T = T - + M 

w ^w 



where: M = ( 4T/Z1Y)^ 4Y 

Constant Wall Temperature with One Wall Insulated . For this case 
the coordinate system was redefined so thb X axis coincided with the lower 
wall. The lower wall boundary condition was taken to be: 



t 



w 



1 



A 

o 




= N 



while zero slope at the upper wall was obtained by equating the wall 
temperature to the calculated temperature one step from the wall: 



Constant Wall Temperaturei(Asymm;etrjo.Ga5e).. With the axis system 
the same as for the insulated wall case, thermal boundary conditions were 
established as; 






Wi 



= N 



1 



tw^Ar. = T = 

^2 o w^ 2 



where and N 2 are positive integers and Nj / N 2 . 
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I 

} 






] 



I 



Noting that: 




^wiAq - 1 _ Nj - 1 






Values of Nj - 1/ N 2 - 1 cf 1/2, 1/3, 1/4 were investigated. 

E . Computational Procedure 

The finite difference solution was obtained by using a CDC 1604 
digital computer to solve the equations at intervals of^X = 0.001 and 
AY =0.1. Velocity distributions were obtained by interpolating plotted 
curves of Bodoia's tabulated results. 

Since the entire solution is dependent upon the first calculations 
at the entrance, starting values were obtained by reducing the grid size 
to/iX = 0.0001 and 4 Y =0.05 for the fipt 20 streamwise stations from 



velocity component (V) was introduced into the energy equation (11) at 
each grid point. The energy solution was then obtained at intervals of 
Y between centerline and the wall at each streamwise station. 

Solutions for each thermal boundary condition were obtained for 
Prandtl numbers of: 0.5, 0.7, 1.0, 1.6, 3.2 and 10. 

Solutions for Prandtl numbers less than 0.5 are excluded by the grid 
size. This may be seen from inspection of the coefficient of the second 
temperature term in (11), which indicates the interdependence of PR,2iX, 



X = 0 to 0.0020. 



The continuity equation (2) was evaluated .artd the cross-wise 
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( 






and AY. In order to obtain convergence, PR must be large enough to 
prevent this term from becoming negative . Thus low Prandtl numbers , 
approaching those of liquid metals, cannot be effectively handled by 
finite-difference methods owing to the large number of calculations 
required by the small grid. 

The computer program, written in Fortran 60, is included in the 
Appendix , 
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Ill, RESULTS 



Tables I through VI contain computed and extrapolated values of 
the heat transfer rates as a function of distance from the entrance for each 
of the thermal boundary conditions for a Prandtl number of 0.7. Figures 
5 through 29 contain cur\'es of the transfer rate versus downstream position 
and temperature and velocity profiles for Prandtl numbers of 0.5, 0.7, 1.0, 
1.6, 3.2 and 10. 

With the exception of the Asym;metric Constant Wall Temperature 
Case, which is expressed in terms of the dimensionless heat flux parameter, 

'ik 

q ; the heat transfer rates are presented as local Nusselt number variation 
with the length parameter, L. Whereas, the computations were carried out 
using the dimensionless length parameter X, which is normally employed 
for hydrodynamic developments; L, based on the hydraulic diameter is 
customarily employed in heat transfer. The conversion is: 

X = 16L 

The use of all finite starting length" in the finite difference equations 
results in a finite initial value for the Nusselt number at X = 0 rather than 
the theoretically infinite value at the entrance. This causes a nearly flat 
slope near the entrance and an artificial inflection point in the Nusselt 
number versus length curve. The inflection point occurs progressively 
farther dowmstream with the increasing Prandtl number, thereby limiting the 
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. 



I 



the useful span of the curves. A reduction of grid size from ZlX = .001 , 

^Y = .lto^^X = .0001, AY - .05 was observed to double the initial 
Nusselt number from 40 to 80, arid move the inflection point upstream. 

This increased the useful span. However, the approach obviously reaches 
a practical limit very rapidly due to the number of calculations and the 
accuracy required of the hydrodynamic solution very close to the entrance . 

In general, for the three thermal boundary conditions in which the 
Nusselt number was used as a heat transfer parameter, the curves for 
Prandtl numbers of 0.5, 0.7, 1»0 and 1„6 behave as expected for values 
of L> lO”"^. For Prandtl numbers of 3.2 and 10, L = 10“^ is the approxi- 
mate lower limit of validity. 

Constant Heat Input . A comparison of local Nusselt numbers obtained 
in the present analysis with those reported by Siegel and Sparrow/^, and 
Han /6/, contained in Fig. 5, indicate close agreement with the results of 
Siegel and Sparrow. Their solution employed a simplified energy equation 
in which the temperature distribution in the boundary layer was expressed 
as a series of polynomials in the transverse coordinate using the down- 
stream station as a parameter. The hydrodynamic islolution employed 
Schiller's approximation. 

The comparison with Han's integro -numerical solution is not as close 
as might appear in Fig. 5 since this curve is for a Prandtl number of 0.8 
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and should lie above the other two curves which are for a Prandtl number 
of 0.7. 



The computed curves of local Nusselt number versus length contained 
in Fig . 6 . for Prandtl numbers of 0.5, 0.7, 1.0, 1.6, 3.2, and 1 0 , i may be 
approximated within 1.0 per cent , for 0.5^ PR ^ 1.6, and L > 10”^ by the 
following; 




8.24 + 



.0186 

1 + .0178 (PR/I.) 



Velocity and temperature profiles for constant heat input are contained 
in Figs . 7 through 10. 

Constant Wall Temperature (Symmetric Case) . Local values of Nusselt 
numbers for this boundary condition are compared in Fig. 11 with those of 
Sparrow /sy, which were computed using the Karma n -P oh la u sen method and 
Schiller's hydrodynamic solution. Again close correlation is indicated 
between the two results. 

The local Nusselt number curves of Fig. 12 may be approximately 
represented withifilO per cent for 0 . 5 PR * 1 . 6 , and L > 10 ^ by the 
following; 



N 



u. 



7.60 + 



.0148 (^Vl) 

1 + .0143 (PR/l)-6 



Velocity and temperature profiles for this boundary condition are 
contained in Figs. 13 through 16. 
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One Wall Constant Temperature - One VV^ll Insulated . The results 
of this computation are presented in Fig. 17. The local Nusselt number may 
be approximated within 10 per cent for 0.5^ PR ^ 1.6, and L > 10~^ by; 



N, = 4.84 + 



0155 



( ) 



1 + .012 

The velocity and temperature distributions are illustrated in Figs. 

18 through 21. 

Constant Wall Temperature (Asymmetric Case) . Dimensionless heat 
flux values for this boundary condition are presented in Figs. 22, 23 and 24 
for Ni - 1 ratios of 1/2, 1/3, and 1/4. 

N2-I 



Note that the upper (hotter) wall curve is smoothly asymptotic to 
the value predicted for the steady state condition. The lower curve 
exhibits a much shallower initial slope, passes through an inflection 
point, and then, approaches a negative asymptotic value of the same 
magnitude and at the same rate as the upper wave , 

Velocity and temperature profiles for Ni - 1 = 1/4 are presented 

N 2 - 1 

in Figs. 25 through 28. 

The behavior of the local Nusselt number for this asymmetrical case 
is depicted in Fig. 29. The lower curve, representing the local Nusselt 
number variation for the lower (cooler) wall, does not depict the actual 
physical situation. 
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Both the denominator and the numerator of; 



N. 



u. 



^ ( Tw * Tw-Ay) 
AY 

T - T 

' w * m. 



( 18 ) 



approach zero at different rates. The curve tends to positive infinity as 
the denominator approaches zero, returns from negative infinity, and 
becomes zero as the numerator goes to zero; and finally approaches the 
actual steady asymptotic value. 



18 



IV. CONCLUSIONS 



From the results of this study it may be concluded that: 

1. The finite difference method, as employed in this analysis, 
yields results which are consistent with previous solutions 
for other geometries /_%J and which are asymptotic to the 
known fully developed values. The method is, however, 
limited to a small range of Prandtl numbers near unity by the 
small grid sizes necessary outside this range. For small 
Prandtl numbers the grid must be very small to obtain con- 
vergence, while for large Prandtl numbers practical grid 
dimensions result in significant inaccuracies near the 
entrance due to the "finite starting length". 

2 . The results of the approximate integral methods suggested 
by Siegel and Sparrow/ 7/ for constant heat input, and by 
Sparrow /5_y^ for constant wall temperature, compare very 
favorably with the finite difference solutions reported here. 
This additional substantiation of the methods suggested by 
these authors is significant because of the practical utility of 
their approach, which requires relatively few calculations, is 
capable of handling a wide range of Prandtl numbers , and can 
predict heat transfer rates very close to the entrance. 
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TABLE I 



CONSTANT HEAT INPUT 
PR = 0 o 7 



X 


L 


Nn 


T 

w 


0016 


oOOOl 


46o2 


lo09 


0032 


.0002 


30.0 


1.14 


0064 


.0004 


23.7 


1.18 


0128 


.0008 


18.0 


1.24 


0160 


.0010 


16.6 


lo26 


0320 


.0020 


12.9 


1.35 


0640 


.0040 


10o3 


1.47 


1280 


.0080 


8.87 


1.62 


1600 


.0100 


8.62 


lo67 




Asymptotic to 


8.24 
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TABLE II 



CONSTANT WALL TEMPERATURE 
(SYMMETRIC CASE) 

PR = 0 . 7 



X 


L 




0016 


.0001 


34.0 


0032 


.0002 


26.0 


0064 


.0004 


19.9 


0128 


- .0008 


15.3 


0160 


.0010 


13.7 


0320 


.0020 


10.5 


0640 


.0040 


8.58 


1280 


.0080 


7.78 


1600 


.0010 


7.63 




Asymptotic to 


7.60 



22 



TABLE III 



ONE WALL CONSTANT TEMPERATURE - ONE WALL INSULATED 

PR = 0 . 7 





L 


Nu 

“x 


T 

w 


0016 


.0001 


33.0 


1.000 


0032 


.0002 


to 

0 

o 


1.000 


0064 


.0004 


19.2 


1.000 


0128 


.0008 


13.0 


1.000 


0160 


.0010 


11.6 


1.000 


0320 


.0020 


9.05 


1.000 


0640 


.0040 


7.20 


1.000 


1280 


.0080 


5.93 


1.001 


1600 


oOlOO 


5.55 


1.021 


4000 


.0250 


4.96 


1.139 




Asymptotic to-' 


4.84 


2.000 
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TABLE IV 



CONSTANT WALL TEMPERATURE 

(ASYMMETRIC CASE) 

- 1 = 1/2 PR = 0 . 7 
N2 - 1 



X 


L 


★ 

^1 


★ 

^2 


.0032 


.0002 


26.2 


57.5 


00064 


.0004 


17.3 


41.4 


.0128 


.0008 


11.9 


28.7 


.0160 


.0010 


10.7 


25.9 


.0320 


.0020 


8.08 


18.7 


.0460 


.0040 


6.18 


13,7 


.1280 


.0080 


4.73 


10.1 


.1600 


.0100 


4.24 


9.15 


.3200 


.0200 


2.24 


6 o 35 


.6400 


.0400 


-.19 


3 o 81 


.9600 


.0600 


- 1.23 


2.76 




Asymptotic to 


- 2.00 


+ 2.00 



24 



TABLE V 



CONSTANT WALL TEMPERATURE 
(ASYMMETRIC CASE) 





- 1 = 1/3 

N2 - 1 


PR = 0.7 




X 


L 


* 

^1 


* 

^2 


0032 


.0002 


28.0 


86.0 


0064 


.0004 


17.3 


62.0 


0128 


.0008 


11.8 


43.0 


0160 


.0010 


10.8 


38.9 


0320 


.0020 


8.09 


28.1 


0640 


.0040 


6.18 


20.6 


L280 


.0080 


4.69 


15.2 


1600 


.0100 


4.15 


13.8 


3200 


.0200 


1.65 


9.84 


6400 


.0400 


-1.58 


6.41 


9600 


.0600 


-2.97 


5.02 




Asymptotic to 


-4.00 


+4.00 
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TABLE VI 



CONSTANT WALL TEMPERATURE 
(ASYMMETRIC CASE) 





N 1 - 1=1/4 
N2 - 1 


PR = 0.7 




X 


L 


* 

^1 


ic 

^2 


.0032 


.0002 


28.2 


116.0 


.0064 


.0004 


16.8 


82.9 


.0128 


.0008 


12.2 


57.7 


.0160 


.0010 


10.5 


50.4 


.0320 


.0020 


8.09 


37.5 


.0640 


.0040 


6.17 


27.5 


.1280 


.0080 


4.66 


20.3 


.1600 


.0100 


4.06 


18,4 


.3200 


.0200 


1.06 


13,3 


.6400 


.0400 


-2.96 


9.03 


.9600 


.0600 


-4.71 


7.27 




Asymptotic to 


-6.00 


+6.00 
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DESIGNATION OF MESH POINTS 
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(b) 

AXIS SYSTEM, SYMMETRICAL CASES', 
CONSTANT WALL TEMPERATURE 
CONSTANT HEAT INPUT 



X 

(C) 

AXIS SYSTEM, NON- SYMMETRICAL CASES', 

ONE WALL insulated; 

CONSTANT, BUT DIFFERENT WALL TEMPERATURES 
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COMPARISON OF VELOCITY DISTRIBUTIONS 
BODOIA AND SCHLICHTING 
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FIGURE 7 

TEMPERATURE AND VELOCITY PROFILES 
CONSTANT HEAT INPUT 
I X=0.005 L=0.0003I25 
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FIGURE 8 



TEMPERATURE AND VELOCITY PROFILES 
X= 0.050 L = 0.003125 
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FIGURE 9 



TEMPERATURE AND VELOCITY PROFILES 
X= 0.250 L =0.0156 
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R6URE 10 

TEMPERATURE AND VELOCITY PROFILES 

X* 1.00 L *0.0625 
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FIGURE 13 

TEMPERATURE AND VELOCITY PROFILES 
I CONSTANT WALL TEMPERATURE 

X=0.005 L =0.0003125 
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FIGURE 14 



TEMPERATURE AND VELOCITY PROFILES 
CONSTANT WALL TEMPERATURE 
X= 0.050 L =0.003125 
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FIGURE 15 

TEMPERATURE AND VELOCITY PROFILES 
, CONSTANT WALL TEMPERATURE 
X= 0.250 L =0.0156 
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. X= 1.00 



FIGURE 16 

AND VELOCITY PROFILES 

WALL TEMPERATURE 
L =0.0625 
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FIGURE 18 

TEMPERATURE AND VELOCITY PROFILES 
ONE WALL CONSTANT TEMPERATURE 
ONE WALL INSULATED 
X= 0.005 L = 0.0003125 
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TEMPERATURE AND VELOCITY PROFILES 
ONE WALL CONSTANT TEMPERATURE 
, ONE WALL INSULATED 

X=0.050 L =0.003125 

) 
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FIGURE 20 

TEMPERATURE AND VELOCITY PROFILES 
ONE WALL CONSTANT TEMPERATURE 

ONE WALL INSULATED 
X= 0.250 L =0.0156 
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FIGURE 21 

TEMPERATURE AND VELOCITY PROFILES 
ONE WALL CONSTANT TEMPERATURE 
ONE WALL INSULATED 
X= LOO L =0.0625 
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FIGURE 25 

TEMPERATURE AND VELOCITY PROFILES 
CONSTANT WALL TEMPERATURE 
(ASYMMETRIC CASE) 

•* X= 0.005 L =0.0003125 

\ 
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FIGURE 26 



TEMPERATURE AND VELOCITY PROFILES 
CONSTANT WALL TEMPERATURE 
(ASYMMETRIC CASE) 

X= 0.050 L =0.003125 
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FIGURE 27 

TEMPERATURE AND VELOCITY PROFILES 
CONSTANT WALL TEMPERATURE 

(ASYMMETRIC CASE) 

X= 0.250 L =0.0156 
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FIGURE 28 

TEMPERATURE AND VELOCITY PROFILES 
CONSTANT WALL TEMPERATURE 
(ASYMMETRIC CASE) 

X= 1.00 L= 0.0625 

1 
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LOCAL NUSSELT NUMBER 
CONSTANT WALL TEMPERATURE 
' (ASYMMETRIC CASE) 

N, -I _ I 
N 2 -I " 4 

PR » 0.7 
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APPENDIX A 

SAMPLE FORTRAN 60 PROGRAM 
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. . J080191F,LUNDBERG 

PROGRAM HEATX23 
CONSTANT HEAT INPUT 
BODOIA VELOCITY DATA 

THIS PROGRAM USES A REDUCED GRID OF. X = .000 1 , Y = . 05 
FOR THE FIRST 20 STREAMWISE CALCULATIONS. 

GRAPH OUTPUTS OF NUSSELT NUMBER VERSUS LENGTH PARAMETER 
FOR'PRANDTL NUMBERS OF 0 . 5 >0 . 7 ♦ 1 , 0 » 1 . 6 , 3 . 2 »AND 10* '■ 

ARE OBTAINED TO 50 PER CENT OF THE HYDRODYNAMIC 
DEVELOPMENT. PRINT OUTPUT INCLUDES NUSSELT NUMBER. 

SOPE OF THE TEMPERATURE PROFILE AT THE WALL . AND THE 
SLOPE OF THE TEMPERATURE PROFILE AT THE WALL, AND THE MEAN 
MIXED MEAN TEMPERATURE. IN ADDITION, VELOCITY AND 
TEMPERATURE VALUES ARE PR.INTED AT EACH GRID POINT AT 
SELECTED STREAMWISE LOCATIONS. 

DIMENSION T (30,23) ,U( 30,23 ) ,V(30, 23) ,Y( 23) ,X( 500) . 

1,XN(500),ITITLE(12) 

READ 500, ( ITITLEI I ) , 1 = 1,6) 

READ 500, ( ITI,TLE( I ) , 1=7,12) 

500 FORMAT (6A8) 

DO 750 KK=1,6 

GO TO(40,A1,42»A3,A4,A5) ,KK 

40 PR=0.7 

GO TO 46 ‘ 

41 PR=0.5 
GO TO 888 

GO TO 46 - ■ . • 

42 PR=1.0 
GO TO 46 

43 PR=1.6 

GO TO 46 , 

44 PR=3.2 
Go TO 46 

45 PR=10. 

46 PRINT 523 

523 FORMAT (IHl) • 

PRINT 17, PR 

17 FORMAT (/////10X,17H PRANDTL NUMBER =,F5.3) 

L = 1 

GO TO 150 

700 READ 701, ((U(I,J), J=13,22)» 1=2,22) 

701 FORMAT (10F6.4) 

DO 710 1=2,22 
DO 711 J=l,12 

711 U( I ,J)=U( I , 13) ' 

710 CONTINUE 
DELX=.0001 
DELY=0.05 
R=1000. 

DO 110 L=2,21 
150 UTM=0.0 
UM=0.0 
XL = L 

X(L)=(XL-1.0)*0.0001 
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DO 111 J=3,22 
XJ = J 

Y(,L) = (XJ-3.0)*0.05 
Y(23)=1.0 
U( 1»J)=1. 

U(l*23)=1.0 

V(1,J-1)=0.0 ■ • 

V( 1 ,J) =0. 

T(1,J)=1.0 

T(1,J+1)=1.0 

T(1,J-1)=1.0 

T(L,1)=T(L»5) 

U(L»1)=U(L»5) 

U(L,23)=0.0 

U(L,2)=U(L»4) 

V(L,3)=0.0 

V(L,2|=V(L»4) 

C THERMAL BOUNDARY CONDITIONS 

T(L»23) =T(L,a2)+.05 
T(L*2)=T(L»4) 

IF(L-l) 109,109»108 

108 V(L»J)= + DELY/(2.*DELX*R) * ( U ( L-1 » J ) -U ( L + 1 ♦ J ) ) ■ 

109 A= DELX/ (U(L,J)*PR*DELY**2) 

V(L,23)=0,0 

B= (DELX*R*V(L, J) ) / { 2 .*DELY*U ( L ♦ J ) ) 

/ . T(L+1,J) = (A-B)*T.(L»J + 1) + (1.-2.*A)*T(L»J) + (A+B)*T(L,J-1,> 

IF(J-3) 113,113,112 

113 UTM=UTM+.5*U(L,J)*T( L»J) 

UM=UM+.5*U( L, J) 

GO TO 114 

112 UTM=UTM+U(L»J)*T(L».J) 

UM=UM+U{L»J) 

114 TM=UTM/UM 
111 CONTINUE 

Tl= ( T(L,23)“T(L,22 ) ) / DELY 
T2= T(L,23) - TM 

XN(L )=4.0*T1/T2 

131 PRINT 118, X(L) »XN(L) ♦T2»T1,TM 

118 FORMAT (// 4H X = ♦ F6 . 4 ♦ 5X ♦ 14H NUSSELT NO. =.F8.5»5X»8H TW-TM =. 
1F8.5,5X»8H SLOPE =»F8.5»5X,5H TM =,F8.5) 

PRINT119, (T(L,J) ♦ J=l»12) 

PRINT119, (T(L»J)» J=13»23) 

119 FORMAT (/9H TEMP. =»12F7.4) 

PRINT120, • (U(L» J) » J=l»12) 

PRINT120, (U(L»J)» J=13»23) 

120 FORMAT (/ 9H X VEL. =,12F7.4) 

PRINT121. (V(LtJ)» J=l»12) 

PRINT121* (V(L,J), J=13»23> 

121 FORMAT (/ 9H Y VELv =»12F7.4) 

IF(L-I) 700,700.110 

110 CONTINUE 

CALL DRAW ( 2 1 »X ,XN ,MOD »0 ,LAB » I T I TLE ♦ .04 , 10 . , 0 , 0 ,0 ♦ 0 , 6 » 6»1»LAST) 

DO 720 N=2,ll 

M=2*N-1 
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/ 



( 



771 

720 

85 



U( 1 .N) =U( 11 »M) 

U(2 »N)=U(21,M) 
T(2,N)=T(21,M) 

PRINT 771, M, U(2,N) ,T(2,N) 
FORMAT ( I10,2F10.-5 ) 

CONTINUE 

L=2 



DO 10 M=2,250 
70 XM=M 



Y=0.0 
DELY=0.1 
DELX=0.001 
X(M)=DELX*XM 
I=L + 1 



N3N=M 
N = L + 1 



82 IF(M- 99) 72,73,73 

72 READ 1 (U( I ,11) ,U( I , 10 ) ,U( I ,9) ,U( I ,8) ,U( I ,7) ,U( I ,6) 
1 ,U( I ,5) ,U( I ,4) ,U( I ,3 ) ,U( I ,2 ) ) 

1 FORMAT (10F6.4) 

73 IFIM-99 ) 50,51,52 
52 IFIM-149) 51,53,54 
54 IF(M-199) 53,55,55 
51 UU=XM- 99. 

U(N,2) =1.4388+.000740*UU 

U(N,3) =1.4292+.000674*UU • • 

•U(N,4) =1.3993+.000492*UU 
U(N,5) =1.3454+.000238*UU 
U(N,6) =1.2628-.000034*UU 
U(N,7) =1. 1467-.000262*UU 
U(N,8) = ,9938-.000410*UU 
U(N,9) = .8022-.000452*UU 



U(N,10)= .5720“.000388*UU 
U(N,11)= .3042-.000230*UU 
GO TO 50 
53 UU=XM-149. 

. U(N,2) =1.4758+.000290#UU 

U(N,3) =1.4629+,000266*UU 
U(N,4) =1.4239+.000194*UU 
U(N,5) =1.3573 + .000092^^UU 
U(N,6) =1.2611-. 000014^fUU 
U(N,7) =1.1336-. 000104*UU 
U(N,8) = .9733-.000160#UU 
U(N,9) = .7796-.000176*UU 
U(N,10)= .5526-.000150*UU 
U(N,11)= .2926-.000092*UU 



GO TO 50 
55 UU=XM-199. 

U(.N,2) =1.4903+. 0000121*UU 
U(N,3) =1.4762+.0000110*UU 
U(N,4) =1.4336+.000b080*UU 
U(N,5) =1.3619+. 0000026^^UU 
U(N,6) =1.2604-. 0000005*UU 
U(N,7) =1.1284-. 0000042*UU 
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U(N,8) = .9653-.0000066*UU 
U(N»9) = .7708-.0000072*UU 
U(.N»10)= .5A51-.0000063*UU 
U(N»11)« .2880-.0000475*UU 
50 UTM=0.0 
UM=0.0 

DO 11 J=2,ll 
XJ=J 

Y( J)«DELYttXJ-.2 
U(L,12)=0.0 
V(L»2)=0.0 
V(L»1)=V( I ,3) 

Y(12)=1.0 

T(L»12)=T(L»11)+0.1 

T(L,1)=T{L,3) 

U(L,1)=U(L»3) 

8 V(L»J)= V(L»J-1) + DELY/(2.*DELX*R) * ( U { L-1 * J ) -U ( L + 1 » J ) ) 

,A= DELX/ (U(LiJ)*PR*DELY**2) 

V(L»12)=0.0 , 

B= (DELX*R«V( L, J) ) / ( 2 .*DELY*U ( L» J ) ) 

T(L+1,J) = (A-B)^^T(L»J + 1 ) + (1.-2.*A)*T(L*J) + <A+B)*T(L,J-1) 

IF(J-2) 13,13»12 

13 UTM=UTM+.5*U(L»J)*T(L»J) 

UM=UM+.5*U( L» J) 

GO TO 14 

12 UTM=UTM+U(L,J)*T(L, J) . . 

UM=UM+U(LtJ) 

14 TM=UTM/UM 
11 CONTINUE 

Tl= ( T(L»12)-T(L»11) ) / DELY 
T2= T(L,12)-TM 

XN(M)=4.0*T1/T2 
IF(N3N-50 ) 31»31»30 

30 IF ( (N3N/25)*25-N3N ) 15»31»31 

31 PRINT 18*X(M) tXNIM) ♦T2»T1.TM 

18 FORMAT {// 4H X = ♦ F6 . 3 t^X ♦ 14H NUSSELT NO. =»F8.5»5X.8H TW-TM =, 

' 1F8.5,5X,8H SLOPE =»F8.5.5X,5H TM =,F8.5) 

PRINT 19, (T(L»J)» J=l»12) 

19 FORMAT (/9H TEMP. =»12F7.4) 

PRINT 20» (U(L,J)» J=l»12) 

20 FORMAT (/ 9H X VEL. =»12F7.4) 

PRINT 21» (V(L,J)» J=1.12) 

21 FORMAT (/ 9H Y VEL. *»12F7.4) 

15 DO 60 J=2,ll 
T(L,J)=T(L+1.J) 

U(L-1»J)=U(L»J) 

60 U(L,J)=U(L+1»J) 

10 CONTINUE ‘ 

IFIKK-6) 760.761.761 
761 MOD=3 
C GO TO 762 

760 M0D=2 

C 762 call’ DRAW(497»X»XN»MOD.O.LAB» ITITLE.. 04.10. ♦0.0»0»0»6» 6»1»LAST). 
750 continue 
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